Primary Analysis
Imputation
Prior to analysis, all outcomes were grouped in 6 batches of similar type to run multiple imputation (MI) on separately, stratified by treatment (Austin 2021). All outcome variables (change scores) in each batch and their respective baseline measures were included in each imputation data set as well as model covariates (e.g. demographics). Additional auxiliary variables were considered to help inform the missing data imputations in each batch (Allison 2009). Candidate auxiliary variables included all baseline, endline, and outcomes (e.g. change scores) across the full dataset as well as a few additional “proxy” variables (e.g. income). Auxiliary variables were included in an imputation data set when two conditions were met: (1) their correlation with the outcomes in that imputation data set were at least |cor| > 0.4 and (2) including those auxiliary variables provided additional observations where the outcomes of interest contained missing values (Allison 2009, Hardt et al 2012, Madley-Dowd et al, 2023). For example, if change scores were missing due to drop out on an outcome variable of interest but baseline data on another variable with complete data were correlated >0.4 (e.g. HEI total score and ADD_SUGARS_ave_bl), this will add precision to the imputation.
The 6 imputation batches/runs and the variables they included are listed below:
Dietary Quality Group 1:
16 Outcomes: Change in ASA24 HEI diet quality scores (total and subscores), Diet ID total score, as well as aMED score.
Covariates: baseline outcome measures, participant’s biological sex at birth, age, race/ethnicity, and education.
Auxiliary variables: ADD_SUGARS_ave_bl
Dietary Quality Group 2:
9 Outcomes: Change in ASA24 average micro and macro nutrients
Covariates: baseline outcome measures, participant’s biological sex at birth, age, race/ethnicity, and education.
Auxiliary variables: HEI2015C13_ADDSUG_bl
Weight Loss:
6 Outcomes: Change in body weight (kg), BMI, percent body weight, as well as three percent body weight achievement outcomes.
Covariates: baseline weight (kg), participant’s biological sex at birth, age, race/ethnicity, and education.
Auxiliary variables: Weightlbs_dietid_change, BMI_dietid_change, BMI_el
Behavioral Group 1:
4 Outcomes: Change in physical activity.
Covariates: baseline outcome measures, participant’s biological sex at birth, age, race/ethnicity, and education.
Behavioral Group 2:
3 Outcomes: Change in self-reported sleep.
Covariates: baseline outcome measures, participant’s biological sex at birth, age, race/ethnicity, and education.
Behavioral Group 3:
15 Outcomes: Change in habit strength (for each behavior assessed, then grouped healthy and unhealthy).
Covariates: baseline outcome measures, participant’s biological sex at birth, age, race/ethnicity, and education.
Multiple imputation was performed in R using the mice
package (van Buuren and Groothuis-Oudshoorn 2011) where the number of
imputations was determined according to the
howManyImputations package (von Hippel 2020). Imputations
were run with the mice() function with 10 iterations each,
using the predictive mean matching (PMM) method.
The line of R code that implements the imputation is:
imputed <- mice(data.imputation, maxit = 10, m = max.runs, predictorMatrix = predM, method = "pmm", print = TRUE)
where max.runs is the largest number of iterations
recommended across all the outcomes of interest in that imputation
run.
Statistical Modeling
For primary outcome of HEI-2015 score:
- Analysis of covariance (ANCOVA) will be used to test baseline to 6-month changes in HEI-2015 total scores for Diet Quality between the WW group and control. Covariates will be included for baseline HEI-2015 total score, participant’s biological sex at birth, age, race/ethnicity, and education.
For secondary outcomes:
For percent weight loss, ANCOVA will be used to test the WW group vs control. Covariates will be included for participant’s biological sex at birth, age, race/ethnicity, education, and baseline weight.
Achievement of 3/5/10 % weight loss is only observed at the 6-month follow-up and is a binary variable. Logistic regression will be used, and covariates will be included for participant’s biological sex at birth, age, race/ethnicity, education, and baseline weight.
The remaining secondary outcomes are measured at baseline and the 6-month post-intervention as continuous scores. ANCOVA will be used to test the baseline to 6-month changes in scores between WW and control groups. Covariates will be included for baseline outcome measures, participant’s biological sex at birth, age, race/ethnicity, and education.
\[ Outcome = Baseline + Sex + Age + Race + Ethnicity + Education + Treatment \]
After analyses were performed for each of the 52 outcomes across 6
imputed data sets, results were pooled across all imputations according
to Rubin’s rules using the pool() function provided in the
mice package in R.
Cohen’s d unadjusted and adjusted effect sizes are
also calculated. Cohen’s d “unadjusted” is a standard measure of d =
(Mean 1 - Mean 2) / (pooled SD) \(=
\text{(Mean 1 - Mean 2)} / \sqrt{\frac{ ( \text{n}_1 - 1) s_1^2 \ + \ (
\text{n}_2 - 1) s_2^2 }{ \text{n}_1 + \text{n}_2 -2 }}\), where
here the mean and SD are calculated on the outcome of interest for each
treatment group, and the d values are averaged across all imputations.
Cohen’s d “adjusted” is based on the model output adjusted for
covariates, and is calculated using the eff_size function
in the R package emmeans (Lenth 2023).
Note: Change scores are calculated as endline minus baseline measurements.
# declare baseline/reference levels for imputed data sets
imputed_long_1 <- imputed_references(imputed_long_1)
imputed_long_2 <- imputed_references(imputed_long_2)
imputed_long_3 <- imputed_references(imputed_long_3)
imputed_long_4 <- imputed_references(imputed_long_4)
imputed_long_5 <- imputed_references(imputed_long_5)
imputed_long_6 <- imputed_references(imputed_long_6)
# make output shell table to save results for each outcome in the analysis loop
cols <- c("Outcome", "Outcome_label", "Mean1", "SE1", "Mean2", "SE2", "Mean_Diff",
"Mean_Diff_LB", "Mean_Diff_UB",
"SE_Diff", "Estimate", "SE", "t", "df", "p_value", "Mean_Skewness", "Min_Skewness",
"Max_Skewness", "Mean_Kurtosis","Min_Kurtosis", "Max_Kurtosis", "Levene_Pvalue",
"Min_Levene","Max_Levene", "Variance_Residual_Ratio", "Variance_Residual_Ratio_Min",
"Variance_Residual_Ratio_Max", "Variance_Outcome_Ratio", "Variance_Outcome_Ratio_Min",
"Variance_Outcome_Ratio_max", "Number of Participants", "Cohens_d_unadjusted", "Cohens_d_adjusted")
output = matrix(nrow = length(outcomes), ncol = length(cols))
colnames(output) = cols
# In order to save OR CI, let's separate ancova and logistic model output:
cols <- c("Outcome", "Outcome_label", "Mean1", "SE1", "Mean2", "SE2", "Mean_Diff",
"Mean_Diff_LB", "Mean_Diff_UB",
"SE_Diff", "Estimate", "SE", "t", "df", "p_value", "Mean_Skewness", "Min_Skewness",
"Max_Skewness", "Mean_Kurtosis","Min_Kurtosis", "Max_Kurtosis", "Levene_Pvalue",
"Min_Levene","Max_Levene", "Variance_Residual_Ratio", "Variance_Residual_Ratio_Min",
"Variance_Residual_Ratio_Max", "Variance_Outcome_Ratio", "Variance_Outcome_Ratio_Min",
"Variance_Outcome_Ratio_max", "Number of Participants", "Cohens_d_unadjusted", "Cohens_d_adjusted")
output.ancova = matrix(nrow = length(outcomes)-3, ncol = length(cols))
colnames(output.ancova) = cols
ancova.outcome_pairs = subset(outcome_pairs, outcomes %in% setdiff(outcome_pairs$outcomes, c("achieve_3_percent_wl", "achieve_5_percent_wl", "achieve_10_percent_wl")))
logistic_cols <- c("Outcome", "Outcome_label", "Prob1", "SE1", "Prob2", "SE2", "OR",
"OR_LCL", "OR_UCL", "OR_SE",
"Estimate", "SE", "z", "df", "p_value", "Mean_Skewness", "Min_Skewness",
"Max_Skewness", "Mean_Kurtosis","Min_Kurtosis", "Max_Kurtosis", "Levene_Pvalue",
"Min_Levene","Max_Levene", "Variance_Residual_Ratio", "Variance_Residual_Ratio_Min",
"Variance_Residual_Ratio_Max", "Variance_Outcome_Ratio", "Variance_Outcome_Ratio_Min",
"Variance_Outcome_Ratio_max", "Number of Participants", "Cohens_d_unadjusted", "Cohens_d_adjusted")
output.logistic = matrix(nrow = 3, ncol = length(logistic_cols))
colnames(output.logistic) = logistic_cols
logistic.outcome_pairs = subset(outcome_pairs, outcomes %in% c("achieve_3_percent_wl", "achieve_5_percent_wl", "achieve_10_percent_wl"))
# originally this was done to save OR CI for the logistic and not ANCOVAs
# but reviewers later asked for Mean diff CIs as well, so could have just done
# one output with same number of columns.
##############
## Modeling ##
##############
# models
# "achieve_3_percent_wl", "achieve_5_percent_wl", "achieve_10_percent_wl", "weight_change"
# get weight_bl as baseline measurement
# Primary will be ITT with multiple imputation
# Formulas could be printed for verification purposes, but it's long.
# save all residuals:
residual_list_dietqual1 <- data.frame()
residual_list_dietqual2 <- data.frame()
residual_list_weightloss <- data.frame()
residual_list_behavioral1 <- data.frame()
residual_list_behavioral2 <- data.frame()
residual_list_behavioral3 <- data.frame()
# ^ used for sensitivity
for (out in outcomes){
# establish which of the imputed data sets contains the outcome of interest:
if(out %in% colnames(imputed_long_1)){
imputed.analysis.data = imputed_long_1
} else if(out %in% colnames(imputed_long_2)){
imputed.analysis.data = imputed_long_2
} else if(out %in% colnames(imputed_long_3)){
imputed.analysis.data = imputed_long_3
} else if(out %in% colnames(imputed_long_4)){
imputed.analysis.data = imputed_long_4
} else if(out %in% colnames(imputed_long_5)){
imputed.analysis.data = imputed_long_5
} else if(out %in% colnames(imputed_long_6)){
imputed.analysis.data = imputed_long_6
}
# setup model formula:
if (out %in% c("changekg_percent_body_wt", "achieve_3_percent_wl", "achieve_5_percent_wl",
"achieve_10_percent_wl" # weight loss
)){ # these outcomes dont match the _change, _bl structure
formula = paste0(" ~ weightkg_bl + Sex_bcf + Age_years + Race2_bcf + Ethnicity_bcf + Education_grouped + Treatment")
} else{formula = paste0("~ ", gsub("_change", "_bl", out), # everything else
" + Sex_bcf + Age_years + Race2_bcf + Ethnicity_bcf + Education_grouped + Treatment")}
# setup model:
if (out %in% c("achieve_3_percent_wl", "achieve_5_percent_wl",
"achieve_10_percent_wl")){
# logistic for binary outcomes
primary.model <- plyr::dlply(imputed.analysis.data, ".imp", function(df)
glm(as.formula(paste(out, formula)), data = df, family = "binomial"))
} else{
# ANCOVA for continuous outcomes
primary.model <- plyr::dlply(imputed.analysis.data, ".imp", function(df)
lm(as.formula(paste(out, formula)), data = df))
}
# Difference between groups:
emm.model <- as.mira(primary.model)
# this command is used so that the emmeans packages recognizes that emm.model is a list of imputed models.
# this is necessary so that the results are pooled using rubin's rules (https://github.com/rvlenth/emmeans/issues/80)
# emmeans is different for logistic vs lm:
if (out %in% c("achieve_3_percent_wl", "achieve_5_percent_wl",
"achieve_10_percent_wl")){
# logistic
emm <- emmeans::emmeans(emm.model, "Treatment", type = "response") # prob, SE, df, asymp.LCL asymp.UCL
# pairs or pairwise
emm.pairs = pairs(emmeans(emm.model, "Treatment", type = "response")) # T1 / T2, OR SE df null z p
# pairs(emmeans(emm.model, "Treatment"), reverse = TRUE) # T2 - T1
# to get OR CI:
emm.pairs.CI = confint( pairs(emmeans(emm.model, "Treatment", type = "response")))
} else{
# ANCOVA
emm <- emmeans::emmeans(emm.model, "Treatment") # estimate, SE, df, z, p
# pairs or pairwise
emm.pairs = pairs(emmeans(emm.model, "Treatment")) # T1 - T2 estiamte SD df z p
# pairs(emmeans(emm.model, "Treatment"), reverse = TRUE) # T2 - T1
# to get OR CI:
emm.pairs.CI = confint( pairs(emmeans(emm.model, "Treatment")))
}
## Diagnostics:
####################
skewness.mean <- mean(moments::skewness(sapply(primary.model, rstandard)))
# calculate skewness
skewness.max <- max( abs( moments::skewness(sapply(primary.model, rstandard)) ) )
skewness.min <- min(moments::skewness(sapply(primary.model, rstandard)))
# save the residuals from the imputed models back onto the data set
# by adding a new column for the standardized residuals from each model
imputed.data.resids <- cbind(imputed.analysis.data,
data.frame(residuals = unlist(llply(primary.model, rstandard), use.names = FALSE)))
# factor treatment to group by later
imputed.data.resids$Treatment = factor(imputed.data.resids$Treatment)
# have to save as several data frames because too large for one
if(out %in% colnames(imputed_long_1)){
residual_list_dietqual1 <- rbind(residual_list_dietqual1,
cbind(imputed.data.resids,
outcome = out))
}
if(out %in% colnames(imputed_long_2)){
residual_list_dietqual2 <- rbind(residual_list_dietqual2,
cbind(imputed.data.resids,
outcome = out))
}
if(out %in% colnames(imputed_long_3)){
residual_list_weightloss <- rbind(residual_list_weightloss,
cbind(imputed.data.resids,
outcome = out))
}
if(out %in% colnames(imputed_long_4)){
residual_list_behavioral1 <- rbind(residual_list_behavioral1,
cbind(imputed.data.resids,
outcome = out))
}
if(out %in% colnames(imputed_long_5)){
residual_list_behavioral2 <- rbind(residual_list_behavioral2,
cbind(imputed.data.resids,
outcome = out))
}
if(out %in% colnames(imputed_long_6)){
residual_list_behavioral3 <- rbind(residual_list_behavioral3,
cbind(imputed.data.resids,
outcome = out))
}
# run the levene's test on each dataset individually, then extract just the p-value
levene <- dlply(imputed.data.resids, ".imp", function(df)
car::leveneTest(residuals ~ Treatment, data = df)$`Pr(>F)`[1])
# mean, min, max
levene.pval <- mean(unlist(levene))
levene.min <- min(unlist(levene))
levene.max <- max(unlist(levene))
# calculate the ratio of residual variance
# first calculate the variance of residuals by group, by imputation.
imputed.data.resids <- imputed.data.resids %>% group_by(.imp, Treatment) %>%
dplyr::summarise(variance.resids = var(residuals), .groups = "keep")
# then group by just imputation to compare the two variances
imputed.data.resids <- imputed.data.resids %>%
group_by(.imp) %>%
dplyr::summarise(ratio.var.resid =
variance.resids[Treatment == rownames(table(imputed.analysis.data$Treatment))[1]]/
variance.resids[Treatment == rownames(table(imputed.analysis.data$Treatment))[2]])
# mean, min, max
variance.residual.ratio.mean <- mean(imputed.data.resids$ratio.var.resid)
variance.residual.ratio.min <- min(imputed.data.resids$ratio.var.resid)
variance.residual.ratio.max <- max(imputed.data.resids$ratio.var.resid)
# and above to take the mean across all imputation sets
# calculate kurtosis mean, max, min
kurtosis.mean <- mean(moments::kurtosis(sapply(primary.model, rstandard)))
kurtosis.max <- max(moments::kurtosis(sapply(primary.model, rstandard)))
kurtosis.min <- min(moments::kurtosis(sapply(primary.model, rstandard)))
# calculate the ratio of outcome variance
# first calculate the variance of outcome by group, by imputation
variance.outcome.ratios <- imputed.analysis.data %>% group_by(.imp, Treatment) %>%
dplyr::summarise(variance.outcome = var(.data[[out]]), .groups = "keep")
# then group by just imputation to compare the two variances
variance.outcome.ratios <- variance.outcome.ratios %>% group_by(.imp) %>%
dplyr::summarise(variance.outcome.ratio =
variance.outcome[Treatment == rownames(table(imputed.analysis.data$Treatment))[1]]/
variance.outcome[Treatment == rownames(table(imputed.analysis.data$Treatment))[2]])
# mean, min, max
variance.outcome.ratio.mean <- mean(variance.outcome.ratios$variance.outcome.ratio)
variance.outcome.ratio.min <- min(variance.outcome.ratios$variance.outcome.ratio)
variance.outcome.ratio.max <- max(variance.outcome.ratios$variance.outcome.ratio)
# Report measure of effect size, Cohen's d
# on imputed for each imputed data set, then average across for all outcomes
# Cohen’s d is calculated from means and SD for differences between groups,
# where d=(Mean1-Mean2)/pooled_SD
# Calculate effect size components by imputation and treatment (for non-binary):
# first calculate components by imputation and treatment
effect_size <- imputed.analysis.data %>% group_by(.imp, Treatment) %>%
dplyr::summarise(Mean = mean(.data[[out]]),
SD = stats:::sd(.data[[out]]),
n = length(.data[[out]]),
.groups = "keep")
# then calculate effect size by imputation (no longer grouping by treatment)
effect_size2 <- effect_size %>% group_by(.imp) %>%
dplyr::summarise(cohensd = (Mean[Treatment == 1] - Mean[Treatment == 2]) /
sqrt( ((n[Treatment == 1] - 1) * SD[Treatment == 1]^2 +
(n[Treatment == 2] - 1) * SD[Treatment == 2]^2) /
(n[Treatment == 1] + n[Treatment == 2] - 2)
),
.groups = "keep")
# Calculate cohens d from model, adjusting for covariates:
effect_size2$cohensd2 <- NA
# emmeans is different for logistic vs lm:
if (out %in% c("achieve_3_percent_wl", "achieve_5_percent_wl",
"achieve_10_percent_wl")){
# for logistic the effect size is the odds ratio
# logistic
# emm <- emmeans::emmeans(emm.model, "Treatment", type = "response") # prob, SE, df, asymp.LCL asymp.UCL
# pairs or pairwise
# emm.pairs = pairs(emmeans(emm.model, "Treatment", type = "response")) # T1 / T2, OR SE df null z p
# pairs(emmeans(emm.model, "Treatment"), reverse = TRUE) # T2 - T1
effect_size2$cohensd2 = data.frame(emm.pairs)$odds.ratio
} else{
for (iter in 1:max(imputed.analysis.data$.imp)){
# ANCOVA
emm_iter <- emmeans::emmeans(primary.model[[iter]], "Treatment") # estimate, SE, df, z, p
# pairs or pairwise
emm.pairs_iter = pairs(emmeans(primary.model[[iter]], "Treatment")) # T1 - T2 estiamte SD df z p
contrast = data.frame(emm.pairs_iter)
effect_size2[iter, "cohensd2"] <- data.frame(eff_size(emm_iter,
sigma = sigma(primary.model[[iter]]),
edf = contrast[, "df"]))$effect.size
}
}
n_participants <- length(rstandard(primary.model[[1]]))
# Save differently for different models:
if (out %in% c("achieve_3_percent_wl", "achieve_5_percent_wl",
"achieve_10_percent_wl")){
# logistic
#emm <- emmeans::emmeans(emm.model, "Treatment", type = "response") # prob, SE, df, asymp.LCL asymp.UCL
# pairs or pairwise
#emm.pairs = pairs(emmeans(emm.model, "Treatment", type = "response")) # T1 / T2, OR SE df null z p
output[which(outcomes == out),] =
c(out, outcome_pairs$outcome_labels[which(outcome_pairs$outcomes == out)],
summary(emm)$prob[1], summary(emm)$SE[1], # prob and SE X
summary(emm)$prob[2], summary(emm)$SE[2], # prob and SE Y
# difference in groups treatment groups:
c(summary(emm.pairs)$odds.ratio, emm.pairs.CI$lower.CL, emm.pairs.CI$upper.CL, summary(emm.pairs)$SE),
unlist(summary(pool(primary.model))[which(summary(pool(primary.model))$term == "Treatment2"),2:6]),
skewness.mean, skewness.min, skewness.max,
kurtosis.mean, kurtosis.min, kurtosis.max,
levene.pval, levene.min, levene.max,
variance.residual.ratio.mean, variance.residual.ratio.min, variance.residual.ratio.max,
variance.outcome.ratio.mean, variance.outcome.ratio.min, variance.outcome.ratio.max,
n_participants,
mean(effect_size2$cohensd), mean(effect_size2$cohensd2) # cohens d
)
output.logistic[which(logistic.outcome_pairs$outcomes == out),] =
c(out, logistic.outcome_pairs$outcome_labels[which(logistic.outcome_pairs$outcomes == out)],
summary(emm)$prob[1], summary(emm)$SE[1], # prob and SE X
summary(emm)$prob[2], summary(emm)$SE[2], # prob and SE Y
# difference in groups treatment groups:
c(summary(emm.pairs)$odds.ratio, emm.pairs.CI$lower.CL, emm.pairs.CI$upper.CL, summary(emm.pairs)$SE),
unlist(summary(pool(primary.model))[which(summary(pool(primary.model))$term == "Treatment2"),2:6]),
skewness.mean, skewness.min, skewness.max,
kurtosis.mean, kurtosis.min, kurtosis.max,
levene.pval, levene.min, levene.max,
variance.residual.ratio.mean, variance.residual.ratio.min, variance.residual.ratio.max,
variance.outcome.ratio.mean, variance.outcome.ratio.min, variance.outcome.ratio.max,
n_participants,
mean(effect_size2$cohensd), mean(effect_size2$cohensd2) # cohens d
)
} else{
# ANCOVA
#emm <- emmeans::emmeans(emm.model, "Treatment") # estimate, SE, df, z, p
# pairs or pairwise
#emm.pairs = pairs(emmeans(emm.model, "Treatment")) # T1 - T2 estiamte SD df z p
output[which(outcomes == out),] =
c(out, outcome_pairs$outcome_labels[which(outcome_pairs$outcomes == out)],
summary(emm)$emmean[1], summary(emm)$SE[1], # mean and SE X
summary(emm)$emmean[2], summary(emm)$SE[2], # mean and SE Y
# difference in groups treatment groups:
c(summary(emm.pairs)$estimate, emm.pairs.CI$lower.CL, emm.pairs.CI$upper.CL, summary(emm.pairs)$SE),
unlist(summary(pool(primary.model))[which(summary(pool(primary.model))$term == "Treatment2"),2:6]),
skewness.mean, skewness.min, skewness.max,
kurtosis.mean, kurtosis.min, kurtosis.max,
levene.pval, levene.min, levene.max,
variance.residual.ratio.mean, variance.residual.ratio.min, variance.residual.ratio.max,
variance.outcome.ratio.mean, variance.outcome.ratio.min, variance.outcome.ratio.max,
n_participants,
mean(effect_size2$cohensd), mean(effect_size2$cohensd2) # cohens d
)
output.ancova[which(ancova.outcome_pairs$outcomes == out),] =
c(out, ancova.outcome_pairs$outcome_labels[which(ancova.outcome_pairs$outcomes == out)],
summary(emm)$emmean[1], summary(emm)$SE[1], # mean and SE X
summary(emm)$emmean[2], summary(emm)$SE[2], # mean and SE Y
# difference in groups treatment groups:
c(summary(emm.pairs)$estimate, emm.pairs.CI$lower.CL, emm.pairs.CI$upper.CL, summary(emm.pairs)$SE),
unlist(summary(pool(primary.model))[which(summary(pool(primary.model))$term == "Treatment2"),2:6]),
skewness.mean, skewness.min, skewness.max,
kurtosis.mean, kurtosis.min, kurtosis.max,
levene.pval, levene.min, levene.max,
variance.residual.ratio.mean, variance.residual.ratio.min, variance.residual.ratio.max,
variance.outcome.ratio.mean, variance.outcome.ratio.min, variance.outcome.ratio.max,
n_participants,
mean(effect_size2$cohensd), mean(effect_size2$cohensd2) # cohens d
)
}
# Optionl returns while loop running:
#print(paste(out, formula))
print(paste0(out, " outcome ", which(outcomes == out), " of ", length(outcomes),"."))
}
# Combined model results
output = as.data.frame(output)
output.ancova = as.data.frame(output.ancova)
output.logistic = as.data.frame(output.logistic)
# Save residuals and output table:
save(residual_list_dietqual1, residual_list_dietqual2,
residual_list_weightloss,
residual_list_behavioral1, residual_list_behavioral2, residual_list_behavioral3,
output,output.ancova, output.logistic,
file = paste0("WINS_Main_Analysis_primary_", Sys.Date(), ".RData"))# load in results from chunk above:
load("WINS_Main_Analysis_primary_2025-01-28.RData")
output = as.data.frame(output)
output.logistic = as.data.frame(output.logistic)
output.ancova = as.data.frame(output.ancova)
# save for wilcox comparison
primary_results <- output
output.ancova.raw <- output.ancova
output.logistic.raw <- output.logistic
# Format columns:
output.ancova[,3:ncol(output.ancova)] = apply(output.ancova[,3:ncol(output.ancova)], 2, as.numeric)
output.logistic[,3:ncol(output.logistic)] = apply(output.logistic[,3:ncol(output.logistic)], 2, as.numeric)
# format p-values:
output.ancova <- output.ancova %>%
mutate(
p_value = case_when(.data[["p_value"]] < 0.001 ~ sub(" ", "", format.pval(.data[["p_value"]], eps = 0.001, digits = 3, nsmall=3)),
TRUE ~ formatC(.data[["p_value"]], digits = 3, format = "f")),
Levene_Pvalue = case_when(.data[["Levene_Pvalue"]] < 0.001 ~ sub(" ", "", format.pval(.data[["Levene_Pvalue"]], eps = 0.001, digits = 3, nsmall=3)),
TRUE ~ formatC(.data[["Levene_Pvalue"]], digits = 3, format = "f")),
Min_Levene = case_when(.data[["Min_Levene"]] < 0.001 ~ sub(" ", "", format.pval(.data[["Min_Levene"]], eps = 0.001, digits = 3, nsmall=3)),
TRUE ~ formatC(.data[["Min_Levene"]], digits = 3, format = "f")),
Max_Levene = case_when(.data[["Max_Levene"]] < 0.001 ~ sub(" ", "", format.pval(.data[["Max_Levene"]], eps = 0.001, digits = 3, nsmall=3)),
TRUE ~ formatC(.data[["Max_Levene"]], digits = 3, format = "f"))
)
output.logistic <- output.logistic %>%
mutate(
p_value = case_when(.data[["p_value"]] < 0.001 ~ sub(" ", "", format.pval(.data[["p_value"]], eps = 0.001, digits = 3, nsmall=3)),
TRUE ~ formatC(.data[["p_value"]], digits = 3, format = "f")),
Levene_Pvalue = case_when(.data[["Levene_Pvalue"]] < 0.001 ~ sub(" ", "", format.pval(.data[["Levene_Pvalue"]], eps = 0.001, digits = 3, nsmall=3)),
TRUE ~ formatC(.data[["Levene_Pvalue"]], digits = 3, format = "f")),
Min_Levene = case_when(.data[["Min_Levene"]] < 0.001 ~ sub(" ", "", format.pval(.data[["Min_Levene"]], eps = 0.001, digits = 3, nsmall=3)),
TRUE ~ formatC(.data[["Min_Levene"]], digits = 3, format = "f")),
Max_Levene = case_when(.data[["Max_Levene"]] < 0.001 ~ sub(" ", "", format.pval(.data[["Max_Levene"]], eps = 0.001, digits = 3, nsmall=3)),
TRUE ~ formatC(.data[["Max_Levene"]], digits = 3, format = "f"))
)
# round non-pvalues to 2 decimal places
p_values = which(colnames(output.ancova) %in% c("p_value", "Levene_Pvalue",
"Min_Levene", "Max_Levene"))
non_pvals <- setdiff(3:ncol(output.ancova), p_values)
output.ancova[,non_pvals] = apply(output.ancova[,non_pvals], 2, function(x) sprintf("%.2f", x))
#round(output.ancova[,non_pvals], digits = 2)
# binary outcomes:
p_values = which(colnames(output.logistic) %in% c("p_value", "Levene_Pvalue",
"Min_Levene", "Max_Levene"))
non_pvals <- setdiff(3:ncol(output.logistic), p_values)
output.logistic[,non_pvals] = apply(output.logistic[,non_pvals], 2, function(x) sprintf("%.2f", x))
rownames(output.ancova) = NULL
rownames(output.logistic) = NULL
# to keep outcome pairs indices, combine CI:
temp = outcome_pairs %>% dplyr::rename(Outcome = outcomes) %>%
dplyr::rename(Outcome_label = outcome_labels)
output.ancova <- merge(temp, output.ancova, by = c("Outcome", "Outcome_label"), sort = FALSE) #all = TRUE,
output.ancova <- output.ancova %>% mutate(
Mean_Diff_CI = paste0("(", Mean_Diff_LB, ", ", Mean_Diff_UB, ")")
) %>%
select(-Mean_Diff_LB, -Mean_Diff_UB) %>%
select(
Outcome, Outcome_label, Mean1, SE1, Mean2, SE2, Mean_Diff,
Mean_Diff_CI, SE_Diff, Estimate, SE, t, df, p_value,
Mean_Skewness, Min_Skewness, Max_Skewness, Mean_Kurtosis, Min_Kurtosis,
Max_Kurtosis, Levene_Pvalue, Min_Levene, Max_Levene,
Variance_Residual_Ratio, Variance_Residual_Ratio_Min,
Variance_Residual_Ratio_Max, Variance_Outcome_Ratio,
Variance_Outcome_Ratio_Min, Variance_Outcome_Ratio_max,
`Number of Participants`, Cohens_d_unadjusted, Cohens_d_adjusted
)Dietary quality
kable(output.ancova[1:14,-c(10, 11, 15:(ncol(output.ancova)-2))],
col.names = c("Variable", "Outcome", "Mean", "SE", "Mean", "SE", "Mean", "95% CI", "SE", "t", "df", "p-value", "unadjusted", "adjusted"),
caption = "Change in ASA24 HEI Diet Quality Scores (Total and Subscores)") %>%
kable_styling(bootstrap_options = "striped", full_width = F) %>% kable_minimal() %>%
add_header_above(c(" " = 2, "Weight Watchers" = 2, "Control" = 2, "Difference"=3, "Model Statistics" = 3, "Cohen's d" = 2))|
Weight Watchers
|
Control
|
Difference
|
Model Statistics
|
Cohen’s d
|
|||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Variable | Outcome | Mean | SE | Mean | SE | Mean | 95% CI | SE | t | df | p-value | unadjusted | adjusted |
| HEI2015_TOTAL_SCORE_change | HEI Total Score | 5.31 | 1.51 | 1.13 | 1.37 | 4.18 | (1.80, 6.56) | 1.21 | -3.45 | 332.56 | <0.001 | 0.30 | 0.37 |
| HEI2015C1_TOTALVEG_change | Total Vegetable | 0.18 | 0.17 | 0.06 | 0.15 | 0.12 | (-0.14, 0.38) | 0.13 | -0.91 | 291.44 | 0.365 | 0.08 | 0.10 |
| HEI2015C2_GREEN_AND_BEAN_change | Greens and Beans | 0.89 | 0.28 | 0.18 | 0.25 | 0.71 | (0.27, 1.15) | 0.22 | -3.19 | 293.87 | 0.002 | 0.21 | 0.35 |
| HEI2015C3_TOTALFRUIT_change | Total Fruit | 0.49 | 0.24 | -0.05 | 0.21 | 0.53 | (0.15, 0.91) | 0.19 | -2.77 | 290.72 | 0.006 | 0.27 | 0.31 |
| HEI2015C4_WHOLEFRUIT_change | Whole Fruit | 0.12 | 0.26 | -0.38 | 0.23 | 0.50 | (0.09, 0.92) | 0.21 | -2.38 | 308.93 | 0.018 | 0.22 | 0.26 |
| HEI2015C5_WHOLEGRAIN_change | Whole Grains | 1.05 | 0.43 | 0.58 | 0.40 | 0.47 | (-0.22, 1.16) | 0.35 | -1.34 | 315.42 | 0.183 | 0.07 | 0.15 |
| HEI2015C6_TOTALDAIRY_change | Total Dairy | -0.52 | 0.39 | -0.24 | 0.35 | -0.28 | (-0.89, 0.34) | 0.31 | 0.88 | 296.63 | 0.381 | -0.02 | -0.10 |
| HEI2015C7_TOTPROT_change | Total Protein Foods | 0.07 | 0.09 | 0.05 | 0.08 | 0.02 | (-0.12, 0.16) | 0.07 | -0.33 | 288.86 | 0.738 | 0.05 | 0.04 |
| HEI2015C8_SEAPLANT_PROT_change | Seafood and Plant Proteins | 0.24 | 0.27 | 0.10 | 0.24 | 0.14 | (-0.28, 0.57) | 0.22 | -0.66 | 315.34 | 0.510 | -0.01 | 0.07 |
| HEI2015C9_FATTYACID_change | Fatty Acids | 1.21 | 0.46 | 0.95 | 0.40 | 0.27 | (-0.44, 0.97) | 0.36 | -0.74 | 282.71 | 0.460 | 0.08 | 0.08 |
| HEI2015C10_SODIUM_change | Sodium | -1.08 | 0.36 | -0.84 | 0.33 | -0.25 | (-0.82, 0.33) | 0.29 | 0.84 | 311.17 | 0.400 | -0.01 | -0.09 |
| HEI2015C11_REFINEDGRAIN_change | Refined Grains | 0.55 | 0.44 | -0.59 | 0.39 | 1.14 | (0.44, 1.84) | 0.36 | -3.19 | 279.74 | 0.002 | 0.33 | 0.36 |
| HEI2015C12_SFAT_change | Saturated Fats | 1.51 | 0.45 | 0.89 | 0.40 | 0.62 | (-0.08, 1.32) | 0.36 | -1.74 | 305.85 | 0.083 | 0.14 | 0.19 |
| HEI2015C13_ADDSUG_change | Added Sugars | 0.60 | 0.28 | 0.44 | 0.25 | 0.16 | (-0.28, 0.60) | 0.22 | -0.72 | 304.15 | 0.470 | 0.02 | 0.08 |
temp <- output.ancova[15,-c(10, 11, 15:(ncol(output.ancova)-2))]
rownames(temp) = NULL
kable(temp,
col.names = c("Variable", "Outcome", "Mean", "SE", "Mean", "SE", "Mean", "95% CI", "SE", "t", "df", "p-value", "unadjusted", "adjusted"),
caption = "Other Dietary Quality Measures") %>%
kable_styling(bootstrap_options = "striped", full_width = F) %>% kable_minimal() %>%
add_header_above(c(" " = 2, "Weight Watchers" = 2, "Control" = 2, "Difference"=3, "Model Statistics" = 3, "Cohen's d" = 2))|
Weight Watchers
|
Control
|
Difference
|
Model Statistics
|
Cohen’s d
|
|||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Variable | Outcome | Mean | SE | Mean | SE | Mean | 95% CI | SE | t | df | p-value | unadjusted | adjusted |
| amed_change | AMED Score | 0.86 | 0.65 | 0.18 | 0.59 | 0.68 | (-0.35, 1.70) | 0.52 | -1.30 | 323.36 | 0.196 | 0.10 | 0.14 |
temp <- output.ancova[16:24,-c(10, 11, 15:(ncol(output.ancova)-2))]
rownames(temp) = NULL
kable(temp,
col.names = c("Variable", "Outcome", "Mean", "SE", "Mean", "SE", "Mean", "95% CI", "SE", "t", "df", "p-value", "unadjusted", "adjusted"),
caption = "Change in Average Micro and Macro Nutrients Between Endline and Baseline") %>%
kable_styling(bootstrap_options = "striped", full_width = F) %>% kable_minimal() %>%
add_header_above(c(" " = 2, "Weight Watchers" = 2, "Control" = 2, "Difference"=3, "Model Statistics" = 3, "Cohen's d" = 2))|
Weight Watchers
|
Control
|
Difference
|
Model Statistics
|
Cohen’s d
|
|||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Variable | Outcome | Mean | SE | Mean | SE | Mean | 95% CI | SE | t | df | p-value | unadjusted | adjusted |
| KCAL_ave_change | Average Total Energy | -468.22 | 70.01 | -241.77 | 63.84 | -226.44 | (-336.93, -115.96) | 56.17 | 4.03 | 336.57 | <0.001 | -0.28 | -0.43 |
| TFAT_ave_change | Average Total Fat | -22.27 | 3.75 | -9.88 | 3.38 | -12.39 | (-18.30, -6.48) | 3.01 | 4.12 | 329.62 | <0.001 | -0.27 | -0.45 |
| CARB_ave_change | Average Total Carbohydrates | -50.49 | 8.50 | -26.79 | 7.72 | -23.70 | (-37.17, -10.23) | 6.85 | 3.46 | 315.60 | <0.001 | -0.27 | -0.38 |
| SODI_ave_change | Average Sodium | -549.96 | 128.91 | -180.12 | 116.19 | -369.84 | (-573.75, -165.92) | 103.63 | 3.57 | 310.54 | <0.001 | -0.27 | -0.39 |
| SFAT_ave_change | Average Saturated Fats | -8.39 | 1.42 | -4.17 | 1.28 | -4.22 | (-6.46, -1.98) | 1.14 | 3.71 | 334.20 | <0.001 | -0.25 | -0.40 |
| SUGR_ave_change | Average Total Sugars | -25.26 | 4.69 | -17.78 | 4.26 | -7.49 | (-14.87, -0.10) | 3.75 | 1.99 | 335.41 | 0.047 | -0.10 | -0.21 |
| ADD_SUGARS_ave_change | Average Added Sugars | -5.01 | 0.97 | -3.18 | 0.87 | -1.84 | (-3.37, -0.30) | 0.78 | 2.36 | 318.86 | 0.019 | -0.15 | -0.26 |
| CHOLE_ave_change | Average Total Cholesterol | -48.79 | 22.80 | -23.50 | 20.25 | -25.29 | (-61.20, 10.63) | 18.25 | 1.39 | 287.58 | 0.167 | -0.07 | -0.15 |
| FIBE_ave_change | Average Fiber | -1.86 | 0.97 | -1.48 | 0.86 | -0.38 | (-1.90, 1.13) | 0.77 | 0.50 | 282.98 | 0.619 | -0.03 | -0.06 |
These measurements were averaged across the ASA24 recalls at endline and baseline, and then the difference was taken (endline - baseline).
Weight loss
temp <- output.ancova[25:27,-c(10, 11, 15:(ncol(output.ancova)-2))]
rownames(temp) = NULL
kable(temp,
col.names = c("Variable", "Outcome", "Mean", "SE", "Mean", "SE", "Mean", "95% CI", "SE", "t", "df", "p-value", "unadjusted", "adjusted"), caption = "Weight Loss Measures") %>%
kable_styling(bootstrap_options = "striped", full_width = F) %>% kable_minimal() %>%
add_header_above(c(" " = 2, "Weight Watchers" = 2, "Control" = 2, "Difference"=3, "Model Statistics" = 3, "Cohen's d" = 2))|
Weight Watchers
|
Control
|
Difference
|
Model Statistics
|
Cohen’s d
|
|||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Variable | Outcome | Mean | SE | Mean | SE | Mean | 95% CI | SE | t | df | p-value | unadjusted | adjusted |
| weightkg_change | Body Weight (kg) | -5.43 | 0.91 | -1.55 | 0.79 | -3.88 | (-5.29, -2.47) | 0.71 | 5.43 | 266.69 | <0.001 | -0.61 | -0.61 |
| BMI_change | BMI | -1.89 | 0.32 | -0.55 | 0.27 | -1.34 | (-1.83, -0.85) | 0.25 | 5.37 | 259.18 | <0.001 | -0.61 | -0.61 |
| changekg_percent_body_wt | Percent Body Weight Change | -5.44 | 0.92 | -1.53 | 0.80 | -3.91 | (-5.35, -2.48) | 0.73 | 5.37 | 257.98 | <0.001 | -0.61 | -0.61 |
temp <- output.logistic[,-c(11, 12, 16:ncol(output.logistic) )]
# combine CIs
temp$CI = paste0("(", formatC(temp$OR_LCL, digits = 2, format = "f"), ", ", formatC(temp$OR_UCL, digits = 2, format = "f"), ")")
temp <- temp[,-which(colnames(temp) %in% c("OR_LCL", "OR_UCL"))]
rownames(temp) = NULL
kable(temp[,c(1:7, 12, 8:11)],
col.names = c("Variable", "Outcome", "Probability", "SE", "Probability", "SE", "OR", "OR CI", "SE", "z", "df", "p-value"), caption = "Logisitc Weight Loss Measures") %>%
kable_styling(bootstrap_options = "striped", full_width = F) %>% kable_minimal() %>%
add_header_above(c(" " = 2, "Weight Watchers" = 2, "Control" = 2, "Odds Ratio" = 3, "Model Statistics" = 3))|
Weight Watchers
|
Control
|
Odds Ratio
|
Model Statistics
|
||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|
| Variable | Outcome | Probability | SE | Probability | SE | OR | OR CI | SE | z | df | p-value |
| achieve_3_percent_wl | Achieved 3% Weight Loss | 0.59 | 0.07 | 0.35 | 0.06 | 2.68 | (1.70, 4.25) | 0.63 | -4.21 | 277.66 | <0.001 |
| achieve_5_percent_wl | Achieved 5% Weight Loss | 0.47 | 0.08 | 0.21 | 0.05 | 3.34 | (2.07, 5.38) | 0.81 | -4.96 | 299.75 | <0.001 |
| achieve_10_percent_wl | Achieved 10% Weight Loss | 0.19 | 0.08 | 0.03 | 0.02 | 7.10 | (3.18, 15.84) | 2.91 | -4.79 | 272.74 | <0.001 |
The effect size from the model, adjusted for covariates, would be the odds ratio of the difference for binary outcomes.
Behavioral
temp <- output.ancova[28:31,-c(10, 11, 15:(ncol(output.ancova)-2))]
rownames(temp) = NULL
kable(temp,
col.names = c("Variable", "Outcome", "Mean", "SE", "Mean", "SE", "Mean", "95% CI", "SE", "t", "df", "p-value", "unadjusted", "adjusted"), caption = "Change in Physical Activity") %>%
kable_styling(bootstrap_options = "striped", full_width = F) %>% kable_minimal() %>%
add_header_above(c(" " = 2, "Weight Watchers" = 2, "Control" = 2, "Difference"=3, "Model Statistics" = 3, "Cohen's d" = 2))|
Weight Watchers
|
Control
|
Difference
|
Model Statistics
|
Cohen’s d
|
|||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Variable | Outcome | Mean | SE | Mean | SE | Mean | 95% CI | SE | t | df | p-value | unadjusted | adjusted |
| METS_change | Total Physical Activity MET | 672.30 | 337.62 | 750.93 | 306.59 | -78.63 | (-600.91, 443.66) | 265.53 | 0.30 | 342.70 | 0.767 | -0.05 | -0.03 |
| sendentary_change | Sedentary | -70.17 | 26.23 | -37.28 | 23.31 | -32.89 | (-74.15, 8.37) | 20.95 | 1.57 | 234.59 | 0.118 | -0.15 | -0.18 |
| moderate_change | Moderate | 231.95 | 211.27 | 293.29 | 191.95 | -61.34 | (-391.88, 269.20) | 168.03 | 0.37 | 330.16 | 0.715 | -0.04 | -0.04 |
| vigorous_change | Vigorous | 511.77 | 227.53 | 515.64 | 205.43 | -3.87 | (-359.30, 351.56) | 180.69 | 0.02 | 332.74 | 0.983 | -0.03 | -0.00 |
temp <- output.ancova[32:34,-c(10, 11, 15:(ncol(output.ancova)-2))]
rownames(temp) = NULL
kable(temp,
col.names = c("Variable", "Outcome", "Mean", "SE", "Mean", "SE", "Mean", "95% CI", "SE", "t", "df", "p-value", "unadjusted", "adjusted"), caption = "Change in Self-Reported Sleep") %>%
kable_styling(bootstrap_options = "striped", full_width = F) %>% kable_minimal() %>%
add_header_above(c(" " = 2, "Weight Watchers" = 2, "Control" = 2, "Difference"=3, "Model Statistics" = 3, "Cohen's d" = 2))|
Weight Watchers
|
Control
|
Difference
|
Model Statistics
|
Cohen’s d
|
|||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Variable | Outcome | Mean | SE | Mean | SE | Mean | 95% CI | SE | t | df | p-value | unadjusted | adjusted |
| sleep_quality_change | Sleep Quality | -0.02 | 0.10 | 0.21 | 0.09 | -0.22 | (-0.38, -0.07) | 0.08 | 2.87 | 309.41 | 0.004 | -0.30 | -0.32 |
| sleep_amount_change | Usual Sleep Amount | -0.08 | 0.05 | 0.02 | 0.04 | -0.10 | (-0.18, -0.02) | 0.04 | 2.48 | 282.76 | 0.014 | -0.23 | -0.28 |
| wake_episodes_change | Wake Episodes | -0.01 | 0.15 | -0.04 | 0.14 | 0.02 | (-0.21, 0.26) | 0.12 | -0.21 | 301.48 | 0.837 | 0.01 | 0.02 |
For sleep quality, lower numbers are better where 1 is very good sleep quality and 5 is very poor sleep quality so higher number of change in sleep quality is worse. Recall for sleep amount 1 is more than usual, 2 is usual, and 3 is much less sleep than usual.
temp <- output.ancova[35:49,-c(10, 11, 15:(ncol(output.ancova)-2))]
rownames(temp) = NULL
kable(temp,
col.names = c("Variable", "Outcome", "Mean", "SE", "Mean", "SE", "Mean", "95% CI", "SE", "t", "df", "p-value", "unadjusted", "adjusted"), caption = "Change in Habit Strength (for each behavior assessed, then grouped healthy and unhealthy)") %>%
kable_styling(bootstrap_options = "striped", full_width = F) %>% kable_minimal() %>%
add_header_above(c(" " = 2, "Weight Watchers" = 2, "Control" = 2, "Difference"=3, "Model Statistics" = 3, "Cohen's d" = 2))|
Weight Watchers
|
Control
|
Difference
|
Model Statistics
|
Cohen’s d
|
|||||||||
|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
| Variable | Outcome | Mean | SE | Mean | SE | Mean | 95% CI | SE | t | df | p-value | unadjusted | adjusted |
| SRBAI_change | SRBAI Habit Strength | 0.87 | 0.31 | 0.41 | 0.15 | 0.46 | (-0.18, 1.11) | 0.31 | -1.48 | 31.15 | 0.149 | 0.48 | 0.51 |
| Avg1_change | Considering Portion Sizes | 1.14 | 0.26 | 0.66 | 0.23 | 0.48 | (0.08, 0.89) | 0.21 | -2.35 | 242.56 | 0.020 | 0.25 | 0.29 |
| Avg2_change | Tracking Food Consumption | 1.46 | 0.29 | 0.91 | 0.26 | 0.54 | (0.08, 1.01) | 0.24 | -2.29 | 248.70 | 0.023 | 0.22 | 0.28 |
| Avg3_change | Consider WW Points | 3.67 | 0.28 | 1.04 | 0.25 | 2.63 | (2.19, 3.08) | 0.23 | -11.67 | 281.27 | <0.001 | 1.28 | 1.35 |
| Avg4_change | Frequency of Eating Vegetables | 0.68 | 0.22 | 0.42 | 0.19 | 0.27 | (-0.07, 0.61) | 0.17 | -1.54 | 263.11 | 0.125 | 0.15 | 0.18 |
| Avg5_change | Frequency of Weighing Self | 1.11 | 0.70 | 0.43 | 0.23 | 0.67 | (-0.72, 2.07) | 0.68 | -0.99 | 26.76 | 0.331 | 0.31 | 0.35 |
| Avg6_change | Frequency of Physical Activity | 0.61 | 0.60 | 0.21 | 0.21 | 0.40 | (-0.82, 1.61) | 0.59 | -0.67 | 28.06 | 0.508 | 0.22 | 0.23 |
| Avg7_change | Talking Kindly to Self After Setback | 0.77 | 0.26 | 0.45 | 0.21 | 0.32 | (-0.11, 0.74) | 0.22 | -1.48 | 203.77 | 0.140 | 0.24 | 0.19 |
| Avg8_change | Arranging Healthy Foods for Easy Access | 1.17 | 0.28 | 0.66 | 0.25 | 0.51 | (0.06, 0.95) | 0.23 | -2.24 | 284.41 | 0.026 | 0.23 | 0.26 |
| Avg9_change | Frequency of Fried Foods | -0.54 | 0.26 | -0.17 | 0.23 | -0.37 | (-0.78, 0.04) | 0.21 | 1.80 | 262.05 | 0.073 | -0.23 | -0.21 |
| Avg10_change | Frequency of Sweets | -1.46 | 0.26 | -0.68 | 0.23 | -0.78 | (-1.18, -0.37) | 0.21 | 3.73 | 281.14 | <0.001 | -0.40 | -0.43 |
| Avg11_change | Frequency of Sugary Beverages | -0.77 | 0.23 | -0.36 | 0.21 | -0.40 | (-0.78, -0.03) | 0.19 | 2.14 | 294.08 | 0.033 | -0.21 | -0.24 |
| Avg12_change | Snacking When Not Hungry | -0.73 | 0.24 | -0.49 | 0.21 | -0.25 | (-0.62, 0.13) | 0.19 | 1.29 | 304.31 | 0.200 | -0.15 | -0.14 |
| UnhSRBAI_change | Unhealthy Grouped | -0.87 | 0.16 | -0.42 | 0.14 | -0.45 | (-0.70, -0.20) | 0.13 | 3.50 | 324.36 | <0.001 | -0.37 | -0.38 |
| healSRBAI_change | Healthy Grouped | 1.31 | 0.13 | 0.58 | 0.12 | 0.73 | (0.53, 0.93) | 0.10 | -7.09 | 311.19 | <0.001 | 0.73 | 0.79 |
